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The ability of high-temperature superconductors (HTSs) to carry very large currents with almost 
no dissipation makes them irreplaceable for high-power applications. The development and further 
improvement of HTS-based cables requires an in-depth understanding of the superconducting vortex 
dynamics in presence of complex pinning landscapes. We present a critical current analysis of a real 
HTS sample in a magnetic field by combining state-of-the-art large-scale Ginzburg-Landau simula¬ 
tions with reconstructive three-dimensional scanning transmission electron microscopy tomography 
of the pinning landscape in Dy-doped YBa 2 Cu 3 07 _ 5 . This methodology provides a unique look 
at the vortex dynamics in the presence of a complex pinning landscape, responsible for the high 
current-carrying capacity characteristic of commercial HTS wires. Our method demonstrates very 
good functional and quantitative agreement of the critical current between simulation and experi¬ 
ment, providing a new predictive tool for HTS wires design. 

PACS numbers: 74.20.De, 74.25.Sv, 74.25.Wx, 05.10.-a 
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1. INTRODUCTION 

Commercial high-temperature superconducting (HTS) 
wires are being successfully applied in a variety of electric 
power equipment for the power grid. Their high current- 
carrying capacity and low dissipation provide multiple 
advantages over conventional conductors.However, 
even commercial second-generation HTS wires have large 
headroom for improvement.^’^ In particular, the critical 
current in these YBa 2 Cu 307_<5 (YBCO) HTS wires de¬ 
creases rapidly in magnetic fields, prohibiting them from 
use in broader applications such as in superconducting 
electrical power generators that could enable lightweight 
and efficient compact systems, e.g., for wind turbines. In 
this work we present the first vortex pinning simulation 
in a reconstructed mixed-pinning landscape obtained by 
three-dimensional (3D) STEM tomography of an actual 
section of a Dy-doped YBCO wire, enabling a new strat¬ 
egy for optimizing the critical current in HTS wires. 

Energy dissipation in superconductors in the presence 
of an applied field arises from the motion of vortices 
driven by the current-induced Lorentz force, ^ thus re¬ 
stricting their mobility through pinning by admixed in¬ 
clusions is the main route to minimize dissipation and 
increase the critical current.At present, the quest 
for higher critical current in HTS is carried out mostly 
via the laborious process of empirical trial and error. 
(Only few systematic studies of the critical current de¬ 
pendence on sizes and densities of defects have been 
published.However, at the fundamental level, the 
basic principles of vortex pinning have been established, 
at least for simple idealized situations, see, e.g., Refs. IS¬ 
IS. A major impediment to rapid progress in improving 


the performance of superconducting wires for applica¬ 
tions is an insufficient understanding of vortex dynamics 
in the complex pinning landscape of real materials. 

Vortex pinning is a complicated collective phenomenon 
controlled by the interaction of vortices with pinning cen¬ 
ters as well as the flexibility of vortex lines and inter¬ 
vortex interactions. The analytical treatment of this 
problem has been limited to qualitative estimates of the 
critical current in simple defect environments such as 
weak pinning by large densities of atomic impurities^^ or 
strong pinning by low densities of strong inclusions. 

On the other hand, numerical simulations of vortex dy¬ 
namics have been used to provide a better insight into 
the pinning mechanisms and to improve the quantita¬ 
tive characterization. Here, the choice among several 
models is determined by a trade-off between complex¬ 
ity and fidelity. Eor example, the minimal approach is 
to consider only the vortex degrees of freedom and treat 
vortices as particles (in thin films) or elastic strings (in 
the bulk). In this case, the dynamics is described by 
an overdamped equation of motion, which takes into ac¬ 
count thermal Langevin forces, see, e.g., Refs. 19-22. 
This approach provides a reasonable description of vor¬ 
tex dynamics at small magnetic fields and small densi¬ 
ties of pinning sites. However, vortex-vortex and vortex- 
defect interactions are treated only approximately in the 
Langevin dynamics approach and vortex cuttings and re¬ 
connections cannot be described within this model. 

In contrast, the time-dependent Ginzburg-Landau 
(TDGL) approach^^’^^ provides an ideal compromise be¬ 
tween approximate phenomenological and exact micro¬ 
scopic descriptions of vortex matter. The TDGL equa¬ 
tions automatically account for the flexibility of the vor- 
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Figure 1: (a) 3D STEM tomogram of a Dy-doped YBCO 

sample (doping level 0.5) processed by the non-uniform illumi¬ 
nation method superimposed with the reconstructed nanopar¬ 
ticles by the IMOD software.(b) Snapshot of the TDGL vor¬ 
tex configuration at an applied magnetic field B = 0.1i7c2 = 
2T and external current J = 0.0052Jdp = 0.2MA/cm^. Iso¬ 
surfaces of the order parameter close to the normal state are 
shown in red and follow both vortex and defect positions. 
The amplitude of the order parameter on the back plane of 
the volume is represented with a blue (normal) to yellow (su¬ 
perconducting) color. 


tex lines, the long-range mutual vortex repulsion, vortex 
cutting and reconnection. Pinning defects of arbitrary 
shape and size can be incorporated via spatial modula¬ 
tion of the critical temperature. Furthermore, the in¬ 
terruption of current paths by replacement of a current- 
carrying superconductor with non-superconducting pin¬ 
ning defects is automatically accounted for by the model. 


TDGL-based numerical simulations have been used sev¬ 
eral times in the past to study various properties of vortex 
matter.^^“^^ 

Recent advancements in computational capabilities in 
combination with efficient parallel solvers for the TDGL 
equation on modern graphics processing units (GPUs)^"^ 
enabled the simulation of rather large 3D samples with 
different pinning landscapes making meaningful predic¬ 
tions for the behavior of critical currents possible. Prac¬ 
tically all theoretical and numerical studies of vortex pin¬ 
ning dealt with idealized models where only one type of 
pinning center is typically considered. However, commer¬ 
cial HTS wires have been engineered with a complex vari¬ 
ety of pinning defects of different sizes and shapes, which 
have empirically been found beneficial for high critical 
current densities. A direct simulation of the vortex be¬ 
havior in such mixed landscape is usually not possible 
because information about internal structures remains 
qualitative and scarce. A recent scanning transmission 
electron microscopy (STEM) tomography study,^^ has 
determined the location and size of inclusions within a 
superconducting Dy-doped YBCO compound used for 
second-generation HTS wires, see Fig. 1(a). This tech¬ 
nique allows for unprecedented exact mapping of inclu¬ 
sion sizes and locations within the sample. 

Here, we combine the advances in large-scale TDGL 
simulations and STEM tomography information to pro¬ 
duce the first comparison of the computed critical cur¬ 
rents for a realistic sample with experimental measure¬ 
ments at different magnitudes and orientations of the ex¬ 
ternal magnetic field. We simulate the complete volume 
and nanoparticles pinning structure of the experimental 
system down to the resolution limit of the tomography. 


2. MODEL 

We use TDGL equations as the main tool for the nu¬ 
merical analysis of vortex dynamics. In the infinite-A 
limit, these equations provide a quantitatively adequate 
description of strong type-H superconductors at high 
magnetic fields. In this limit, the Maxwell-type equation 
for the vector potential is eliminated and we concentrate 
on the remaining equation for the superconducting order 
parameter ^ = '0(r, t), 

{dt + = e{Y)ip - + (V - iAfxl) + C(r, t), (1) 

where /i = /i(r, t) is the chemical potential, A is the 
vector potential associated with the external magnetic 
field B as B = V X A, and is the temperature- 

dependent (i-correlated Langevin term.^^ In Eq. (1), writ¬ 
ten in the dimensionless form, the unit of length is given 
by the superconducting coherence length the unit of 
time is to = 47rcrA^/c^, where A is the London penetra¬ 
tion depth, (7 the normal state conductance, and the unit 
of magnetic field is given by the c-axis upper critical field 
Hc 2 = /ic/2e^^. Here —e is the electron’s charge and c 



























3 


the speed of light. In our notations of the TDGL equa¬ 
tions, the a 6 -plane of the HTS is in the x^-plane of the 
simulation and the c-axis is along the ^-direction. 

To account for the anisotropy and layered structure of 
the material, we replace the 2 ;-component of the Lapla- 
cian term in Eq. (1), (V^ — iAz)‘^ 2 p, with the discrete 
second derivative [ip{z + — 

2^p{z)]/^‘^h‘l^ where 7 is the anisotropy factor and hz is 
the grid point spacing in 2 ;-direction, hz = LzjNz. While 
similar numerical discretizations are used for the in-plane 
directions {x and y), the discreteness in ^ direction is an 
essential part of the model. This description is known 
as the Lawrence-Doniach model.In this model the ef¬ 
fective interlayer coupling is controlled by the parame¬ 
ter ^/^hz- In particular, the layered structure leads to 
intrinsic pinning of vortices for magnetic fields applied 
in-plane. 

The dimensionless function e(r) oc Tc(r) — T vanishes 
at the local critical temperature T ^c(r). The ex¬ 
plicit dependence of the critical temperature Tc(r) on 
spatial coordinate r is used to model large-scale inho¬ 
mogeneities in the superconductor and is a convenient 
way to introduce pinning effects. We use e(r) = 1 inside 
the superconductor. This corresponds to setting the co¬ 
herence length ^ at a given temperature as the unit of 
length. Inside non-superconducting inclusions we choose 
e(r) = -1. 

The total (normal and superconducting) in-plane re¬ 
duced current density in units of Jo = hc^/S7reX‘^^ is 
given by the expression 

J = Im [V^* (V - iA)^] - V/i. ( 2 ) 

The maximum theoretical depairing current density is 
then Jdp = 2/(3\/3)Jo. In addition to Eq. ( 1 ), we solve 
the Poisson equation, VJ = 0, for the scalar potential ji. 

To determine the critical current value, Jc, we apply an 
external current, J, in the a 6 -plane (in the x-direction) 
and ramp it down from the resistive to the superconduct¬ 
ing state of the sample. Each time step we calculate the 
electric field in direction of the applied current averaged 
across the sample, i.e., E = (V^^/i). The I-V curve is 
then obtained by averaging the electric field E over the 
steady states reached after the transient regime follow¬ 
ing each current ramping event. We use a finite electric 
field criterion E^ = 10~^Eq to determine the critical cur¬ 
rent, where Eq = Jq/ a is the electric field unit. In other 
words, we compare electrical field E induced by the ex¬ 
ternal current J with a certain electric field level Eq. 
The defined level is sufficiently low for a reasonable defi¬ 
nition of the critical current. However, it is much higher 
than the level of dissipation corresponding to the value 
of 1 /iV/cm, routinely used as a practical criterion for Jq. 
Therefore, the simulated critical currents are expected to 
be somewhat higher than the experimental ones. 

To simulate the anisotropy of the critical current we 
apply and rotate the magnetic field B from the a 6 -plane 
to the c-axis of the HTS keeping it perpendicular to the 
external current. In the coordinate system of the simula¬ 


tion, the applied field B = H[0, sin 6 >, cosO] can be de¬ 
scribed by the gauge A = xB [O, cos 6 >, — sin 6 >], where B 
is the absolute value of the field and 0 the angle with 
respect to the c-axis. 

In the simulated volume, pinning centers are posi¬ 
tioned in strict accordance with their positions and sizes 
as measured by 3D STEM tomography^^ in the experi¬ 
mental sample. The original reconstruction of the super¬ 
conductor volume together with the particles obtained 
by the segmentation of the STEM tomogram using the 
IMOD software is shown in Eig. 1(a). The defect in¬ 
clusions are confined within a rectangular box of size 
534 X 524 X 129 nm^. This box contains 71 almost spher¬ 
ical particles with sizes ranging from 12.2 to 100 nm. 
Eor the numerical analysis we use a coherence length 
^ = 4.2 nm in the a 6 -plane, which is close to the ex¬ 
perimental value at 77 K and an anisotropy factor of 
7 = 5 suitable for YBCO. We simulate a volume of 
size Lx X Ly X Lz = 128^ x 128^ x 32^, which cor¬ 
responds to 538 X 538 x 134 nm^ with our choice of 
The non-superconducting particles with diameters rang¬ 
ing from 2.90 to 23.8^ occupy about 8.1% of the sim¬ 
ulated volume. A typical vortex configuration for mag¬ 
netic field B = 0 . 1 i 7 c 2 = 2 T is shown in Eig. 1(b), where 
Hc 2 = 20 T is the c-axis upper critical field at 77 K. Iso¬ 
surfaces = 0.1 of the order parameter, shown in red 
in Eig. 1(b), reveal the vortex positions and the contours 
of the pinning landscape. On the back and sides of the 
simulated volume, a color-code indicates the amplitude 
of 


3. RESULTS 

In Eig. 2, we present dependence of the critical current, 
Jc(H), on the magnetic field applied parallel to the c- 
axis of the HTS ( 6 > = 0°). Experimentally,^^ at 77 K and 
in a magnetic field ranging from 0.02 to 1.5 T the crit¬ 
ical current follows a power-law dependence Jc oc B~^ 
with an exponent decreasing with increasing Dy doping 
a ~ 0.74 (green curve) and 0.69 (yellow curve) for Dy 
doping 0.5 and 0.75 respectively.^ At lower magnetic 
field, B < 0.01 T, the critical current is more or less 
independent of the external field, since self-field effects 
dominate. In this work we concentrate on higher-field 
region, B > 0.005iJc2 = 0.01 T, where this effect is neg¬ 
ligible. The red curve in Eig. 2 presents the simulated 
dependence Jc{B) from the 0.5 Dy doped sample shown 
in Eig. 1(a), using only the STEM-resolved nanoparti¬ 
cles. Eor B ranging from 0.005i7c2 (0.1 T) to 0.075i7c2 
(1.5 T), corresponding to the range used in the experi¬ 
mental measurements, the simulated Jc{B) is also de- 


^ Note, that we refitted the experimental data in the same field 
range as for the simulation data. Therefore the exponents are 
slightly different from the ones reported in Ref. 36. 
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scribed by a power law but with the slightly larger ex¬ 
ponent, Of ^ 0.80. To show the effects of defects smaller 
than the resolution of the tomography, we also added 
1500 spherical particles of size 2^ (8.4 nm) to the recon¬ 
structed defects. These background inclusions reduce the 
exponent of the Jc{B) dependence to a 0.68, see blue 
curve in Fig. 2. 

In Fig. 3, we present the anisotropy of the critical 
current in a tilted magnetic field, always perpendicu¬ 
lar to the applied current. The dependencies Jc{0) on 
the tilt angle of the magnetic field 0 with respect to the 
c-axis for three different magnetic fields B are shown. 
The red curves correspond to the simulation with only 
the defects detected by STEM tomography, while the 
blue lines correspond to the simulation with the addi¬ 
tional background inclusions. The anisotropy of Jc{0) 
increases with magnetic field and at the highest fields, 
B > 0.01i7c2, Jc is nearly fiat at all angles except for 
a narrow peak around 0 = 90° caused by the simu¬ 
lated intrinsic pinning due to the layered structure of 
the material. The ratio Jc(90°)/Jc(0°) varies depending 
on the field: Jc(90°)/Jc(0°) 2.9 for B = 0.005iJc2 

[Fig. 3(a)], 3.9 for B = 0.05iJc2 [Fig. 3(b)], and 5.3 

for B = 0.1i^c2 [Fig. 3(c)]. This ratio is mainly deter¬ 
mined by the Lawrence-Doniach parameter ^/^hz = 0.4. 


4. DISCUSSIONS 

The almost quantitative agreement between the Jc{B) 
exponents in the simulation (0.80) and the measurements 
(0.74) is quite striking for this first of its kind numerical 
simulation of a real 3D sample. The small discrepancy 
in the exponents is most likely related to the presence 
of background inclusions that are smaller than the reso¬ 
lution of the STEM tomogram. In fact, we have shown 
that those additional small defects reduce the field ex¬ 
ponent in the simulations. Such small inclusions, that 
several groups are trying to control through chemical 
processes,become especially relevant at high fields 
when all strong pinning sites are already occupied. In 
this regime it was also found in experiments that their 
effect is to reduce the exponent of the field dependence 
of the critical current.Specifically, the addition of 
1500 spherical nanoparticles with diameter 2^, occupy¬ 
ing only 1.2% of the sample volume, reduces the expo¬ 
nent from a ~ 0.80 (red curve in Eig. 2) to a ~ 0.68 (blue 
curve). The latter essentially coincides with the experi¬ 
mental value a ^ 0.69 for 0.75 Dy doping. The increased 
Dy content is expected to create more nanoparticles and 
thus to reduce the exponent a. Our simulation therefore 
appears to correctly reproduce this trend. 

Next, we compare the absolute values of the critical 
current. In the experiment,the critical current density 
varies from Jc = 2.5 MA/cm^ at magnetic field B < O.I T 
to Jc = 0.4MA/cm^ at 5 = I T, for B applied along 
the c-axis. Erom our numerical simulations for these two 
fields we obtain Jc = O.IOJdp = 3.8MA/cm^ and Jc = 


B[T] 
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Figure 2: Critical current Jc as a function of the magnetic 
field B applied along the c-axis of YBCO. The simulated field 
dependence (red curve) with only the nanoparticles observed 
by STEM tomography in the sample with 0.5 Dy doping ex¬ 
hibits almost the same exponent as the experiment (green 
curve). Adding small inclusions, 2^ in diameter, to the sim¬ 
ulation makes the dependence less steep (blue curve), which 
yields an exponent very similar to the experimental one in the 
sample with 0.75 Dy doping (yellow curve). The experimen¬ 
tal Jc{B) dependencies were taken from Refs. 23,36 for B 
ranging from O.IT (0.005iJc2) to 1.5 T (0.075iJc2). Inset: 
Double-logarithmic scale to show the exponents. 

0.021 Jdp = 0.77MA/cm^, respectively. Here we estimate 
^dp = 37MA/cm^ using ^ = 4.2 nm and A = 270 nm at 
77 K. The absolute values are also in a reasonable agree¬ 
ment (within a factor of two) with the experiment. This 
semi-quantitative agreement is quite remarkable consid¬ 
ering that we included only the Dy nanoparticles as pin¬ 
ning centers but neglected other types of defects present 
in YBCO films such as twin boundaries, atomic point 
defects, stacking faults, and dislocations. In addition, 
at least some of the larger simulated values of Jc may 
be attributed to the higher critical electric field criterion 
used for their determination, as explained in the model 
section. This suggests that the nanoparticle defects play 
the dominant role in the field dependence of the critical 
current in the measured commercial tape. 

As shown in Eig. 3(b), the overall shape of the angular 
dependence appears consistent with the experiment,in 
particular the fiat plateau around 0 = 0° (H||c) is clearly 
reproduced. The overall shape of angular dependences is 
expected. As the nanoparticles have an isotropic shape 
and no c-axis correlated pinning centers (e.g. disloca¬ 
tions or twin boundaries) were introduced, there is no 
sharp peak in Jc for B\\c. The only correlated pinning 
is the intrinsic pinning in the a6-plane which does yield 
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Figure 3: Simulated critical current Jc as a function of the angle 0 of the applied magnetic field and for different values of 
magnetic field: (a) B = 0.005i/c2, (b) B = 0.05Hc2 with experimental data^® for 0.5 Dy doping and B = IT superimposed, 
and (c) B = 0.1i/c2- 


a peak in Jq. In addition, in an anisotropic material 
Jc is expected to have smooth angular dependence with 
the typical angle 0 ^ arctany, as we observe at small 
field, Fig. 3(a). With increasing magnetic field, the sim¬ 
ulated Jc anisotropy is found to increase, a trend also in 
agreement with experiments.^ In absolute terms, the sim¬ 
ulated anisotropy is Jc(90°)/Jc(0°) ^ 3.9 for B = IT, 
whereas the experimental anisotropy is somewhat smaller 

2.2 for 0.5 Dy doping and ~ 4.5 for the undoped sam¬ 
ple. This discrepancy is easily explained since our choice 
of the Lawrence-Doniach parameter =0.4 yields 

qualitatively reasonable Jc anisotropy but would need to 
be fine-tuned against experimental measurements of the 
intrinsic pinning alone in clean single crystal for really 
quantitative agreement. More importantly, beyond in¬ 
trinsic pinning we did not take into account other types 
of a6-plane pinning centers, such as flat precipitates, as 
they do not appear in the tomogram. Nonetheless, our 
model of a6-plane pinning is reasonable on a qualitative 
level, and the simulated angular dependence does exhibit 
the proper trends for the combined effects of nanoparti¬ 
cles and intrinsic pinning. 

Finally, from a practical point of view, the most rel¬ 
evant question is probably how “optimal” the observed 
pinning configuration is in terms of absolute critical cur¬ 
rent values. The question can be answered in part by 
looking at a recent study,^^ by some of the present au¬ 
thors, on monodisperse spherical defects using the same 
approach. It was found that the optimal critical cur¬ 
rent for such defects, in magnetic fields ranging from 
0.05i^c2 to 0 . 1 i^c 2 , was achieved for defect diameters 
ranging between 3.5^ to 4^ and occupying a volume frac¬ 
tion of about 20%. In the sample measured by tomog¬ 
raphy, the nanoparticles occupy a total volume fraction 
of 8.1%. The distribution of diameters of these particles 
peaks in the range 3^ to 7^ (72% of all particle are in this 
range) [see histogram Fig. 5(d) in Ref. 23]. However, the 
volume fraction occupied by the above mentioned “op¬ 
timal” particles is only 0.8%. Even counting all parti¬ 
cles from 2^ to 7^ and taking into account the added 
background inclusions, barely increases the volume frac¬ 


tion of these pinning centers to 2%, which is ten times 
less than the optimal density found in Ref. 40. Most of 
the defect volume fraction is actually due to a few large 
nanoparticles, which may capture several vortices at a 
time [see Fig. 1(b)]. According to our simulations, these 
large defects contribute little to pinning and also reduce 
the effective cross-section of the sample, thus negatively 
affecting the critical current. 

In fact, the critical currents found in this study are 
about three times smaller than for the optimal pinning 
configuration for monodisperse spherical defects men¬ 
tioned above. This therefore suggests that the present 
configuration is far from optimal and could be im¬ 
proved by skewing the distribution of nanoparticles to¬ 
ward smaller sizes, for instance 10-20 nm in diameter, 
and simultaneously by raising their density. This is in 
line with recent experimental results that find extremely 
effective pinning with nanoparticles of size 15-30 nm^^ 
or 20-80 nm^^ in YBCO, and 8nm in pnictides.^^ Also, 
similar results were observed with smaller irradiation- 
induced defects, 5 nm in size.^’^^ 


5. CONCLUSIONS 

In summary, we numerically simulated the supercon¬ 
ducting vortex dynamics in a real pinning landscape of 
nanoparticles inside a commercial Dy-doped YBCO tape, 
using the same sample size as in the experiment. The po¬ 
sitions and sizes of these nanoparticles were obtained di¬ 
rectly through STEM tomography reconstruction of the 
sample nanostructure. We obtained good qualitative and 
almost quantitative agreement in the functional depen¬ 
dencies and absolute values of the critical current be¬ 
tween simulation and experiment. One can expect that 
more detailed STEM tomography studies, and in partic¬ 
ular more quantitative characterizations of pinning cen¬ 
ters, should improve the quantitative agreement. In ad¬ 
dition, there is still room for improving the simulation 
in terms of exact treatment of intrinsic pinning as well 
as lowering the voltage criterion towards the one used 
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in the industry. Our results, however, show a promis¬ 
ing pathway for the quantitative analysis and optimiza¬ 
tion of vortex dynamics in various realistic defects envi¬ 
ronment and physical conditions, using a combination of 
3D tomography analysis and large-scale time-dependent 
Ginzburg-Landau simulations. 
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